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ABSTRACT 

Both observations and numerical simulations show that stellar convective motions are composed of 
semi-regular flows of convective rolling cells and the fully developed turbulence. Although the convec- 
tive rolling cells are crucial for the properties of the stellar convection that transports heat and mixes 
materials in the stellar interior, their contributions have not been included in turbulent convection 
models proposed up to now. We simplify the structure of the convective rolling cells as a cellular 
pattern moving circle by circle with different angular velocities around the center, estimating their 
typical size by solving approximately for the temperature difference over the stationary temperature 
background and their average shear of velocity by evaluating approximately their kinetic energy trans- 
formed by themselves working as thermal engines from the heat involved in the convective rolling cells. 
We obtain a steady state solution in the fully local equilibrium which is similar to what is obtained 
in the standard mixing-length theory, by applying such model assumptions to the standard k-e model 
and properly choosing the model parameter c^s. Accordingly, we propose a k-uj model to include the 
transport effect of turbulence in stars. Preliminary results of their applications to the sun and other 
stars with different masses and in different evolutionary stages show good agreements with results of 
the standard mixing-length theory and results of numerical simulations for the stellar convection. 
Subject headings: stars: interiors — stars: evolution — convection — turbulence 



1. INTRODUCTION 

Thermal convection is a common phenomenon in the 
stellar interior. Convective flows are driven in stars 
by the thermal buoyancy, which results from unstable 
density stratifications. Both observations and numeri- 
cal simulations show that convective motions are char- 
acterized by deterministic structures of different scales 
such as unsteady convective rolls and semi-regular con- 
vective cells created by large-scale instabilities of the 
buoyancy, along with the fully developed turbulence 
due to extremely high Reynolds numbers in the stel- 
lar convection zones. The bulk of a convection zone 
is stochastically filled up with numerous unsteady con- 
vective roll s pushing each other and leaving little space 
in be t ween (IMea kin fc Arnettll2067tlKenieres fc Hanialid 
120091: iTramped ach 20ll). Close to the boundaries of 
the convective region, semi-regular convective cells jostle 
onto the surface and roughly line up se en for examples 
as th e Raylcigh-Bcnard convection (Keni eres fc Hanialic 
20001 ^2006 ) and the solar granulation (|Stein fc Nordlund 



1998; Frev tag et al.llMl . 

Such semi-deterministic and well-organized structures 
have already been recognized to have fundamental signif- 
icance on the overall properties of the convection. From 
the dynamical point of view, the buoyancy drives the 
fluids continuously moving up and down alternatively in 
the convective rolling cells by doing the mechanical works 
respectively on the corresponding flows. Due to the equa- 
tion of mass conservation, the flows form a cellular pat- 
tern with the warm fluids moving up in the center and 
the cold fluids moving down along the borders. Such a 
cellular motion creates shears of the velocities between 
the neighbouring layers in the convective rolling cells, 
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which contributes to the generation of turbulence. From 
the thermodynamic point of view, the cellular motion in- 
troduces inhomogeneous temperature perturbations over 
the stationary temperature background in the convec- 
tive rolling cells, which results in extra fluxes transfer- 
ring heat into or out of the boundaries of the convective 
cells. Therefore, each of the convective rolling cells runs 
as an individual heat engine, within which fluids absorb- 
ing heat at the higher temperature boundary while re- 
leasing heat at the lower temperature boundary, and the 
convective cell operates cycle by cycle to transform heat 
into the mechanical work, which is the energy supply of 
the cellular motions themselves in the convective rolling 
cells. 

Thermal convection is an essential input physics for 
the stellar modelling. It not only carries heat outward 
from the stellar interior, but mixes different materials 
completely in the convective regions as well. In order 
to incorporate these effects of the convection in stellar 
models, various convection models have been proposed, 
including the most widely used mixing-length theory 
([Bohm-Vitens&.1953l 119581) and more rec ently turb ulent 
convection models (lXiondll980ll l989: D eng et all 120061: 
Canutolll997t iCanutq fc Dubovikovi, 19981 : iCanuto et al.l 



20011 : iLifc Yang| |2001[). Turbulent convection models are 



based on the averaged hydrodynamic moment equations, 
and include many properties of turbulence such as the 
local effects of generation, dissipation, anisotropy, and 
the nonlocal effect of turbulent transport. However, all 
of them suffer from a common problem, i.e., the contri- 
butions from the convective rolling cells have not been 
considered up to now. Consequently, these models either 
do not result in a steady state solution in the fully local 
equilibrium, or have to adopt a length model of turbu- 
lence similar to what has been introduced in the standard 
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mixing-length theory. 

In order to incorporate properly the effeets of the con- 
vective rolling cells into the turbulent convection models, 
we have tried to find out the qualitative properties of the 
convective rolling cells in the stellar convection zone, by 
approximating the configuration of the convective rolling 
cells based on the results of the numerical simulations as 
a simple eddy structure with fluids moving circle by cir- 
cle around the center of the eddy. The basic equations 
governing the motions of the convective rolling cells are 
given in Sect. 2. In Sect. 3, equations of the standard 
k-e model for turbulence are introduced firstly, and then 
their steady state solutions in the fully local equilibrium 
are discussed. Based on the advanced turbulence mod- 
els, we introduce new models for the rate of buoyancy 
production and for the turbulent heat flux with modifl- 
cations suitable for the turbulent motions in the stellar 
convective regions. With the aid of above models, the 
temperature gradient in the convective region can be de- 
rived. In Sect. 4, we first estimate the typical size of the 
convective rolling cells by solving approximately for the 
temperature difference in a convective eddy, and then 
evaluate the shear production rate of turbulence by a 
model of the velocity shear. For vertical flows similar to 
those down-drafts seen in the numerical simulations of 
the stellar convection, we evahiate their development by 
introducing a decaying distance of the temperature dif- 
ference along the vertical flow. In Sect. 5, we introduce 
a new choice of parameter in the standard k-e model, 
which leads to the well known properties of turbulence 
in the stellar convection zones. By comparing with the 
length model of turbulence in the standard mixing- length 
theory, we introduce a new macro-length model for the 
turbulent convection models. Based on above model as- 
sumptions, we obtain a steady state solution of the k-s 
model in the fully local equilibrium, and compare it with 
the result of the mixing-length theory. In accordance 
with all model assumptions introduced above, we sug- 
gest a k-co model for the stellar turbulent convection in 
Sect. 6, and apply it to the sun and some other stellar 
modc^ls with different masses in Sect. 7. We summarize 
in Sect. 8 our main conclusions, and discuss some im- 
portant effects that do not yet properly incorporated in 
our stellar turbulent convection model. 



2. EQUATIONS FOR MOTIONS OF CONVECTIVE 
ROLLING CELLS 

We approximate the stellar convection zone as a plane 

parallel structure, and establish Cartesian coordinates 
(x, y, z) with the vertical direction z along the stellar 
radius. 

We suppose that the stellar turbulent convection is 
composed of a more or less regular motion of convec- 
tive rolling cells and a fully developed turbulent motion. 
The semi-regular motions of the convective rolling cells 
are determined by the Reynolds-averaged Navier-Stokes 
equation: 



respectively as: 



dt ^ dxi 



1 dp 
p dxi 



In Eq. (1), the velocity of the convective flow Ui and the 
velocity fluctuation of turbulence Ui can be expressed 



Ui = {U,V,W), 
Ui = {u,v,w), 



(2) 



rji the component of gravitational acceleration in the Xi 
direction, p the pressure, p the density, and v the molec- 
ular viscosity. Note that summation should be made over 
all three components if the subscript i appears twice in 
the same expression. 

The stellar convection zone is as a whole in the static 
state, satisfying the equation of hydrostatic equilibrium: 



dp 
'd'z 



= -P09, 



(3) 



where po is the stationary density distribution in the con- 
vection zone and the gravitational acceleration g is de- 
flned as: 

gi = (0, 0, -g) . (4) 

The unstable stratification results in a density differ- 
ence Ap in the fluid, which is related to the temperature 
difference AT under the Boussinesq approximation: 



Ap 

Po 



-I3AT, 



(5) 



where T is the temperature and the thermodynamic co- 
efficient (3 is defined by: 



^ P\dT)^ 



(6) 



Substituting Eqs. (3) and (5) into Eq. (1), the convective 
motion in the steady state is described by: 



dxj 



-09iAT ■ 



duiUj 

dxi 



(7) 



It can be noticed that the first term on the right hand 
side of Eq. (7) is the buoyancy, and the second term 
describes the resistance due to the Reynolds stress. 

The Reynolds-averaged equation of energy conserva- 
tion for the convective rolling cells can be expressed as: 



pT {— + U - — f A— 

dt ^ dxi J dxi V dxi 



(8) 



where s is the entropy of the stellar matter and Fq the 
convective heat flux in the Xi direction. The thermal 
conductivity of radiation A is defined by: 



A = 



3pK ' 



(9) 



where k is the opacity of the stellar matter and a the 
Stefan-Boltzmann constant. 

Usually the total energy fiux F is assumed to be a 
constant in the stellar convective envelope, which results 



F = Fc-X 



dTp 
dz ' 



(10) 



where Tq is the stationary temperature distribution in 
the convection zone and Fc the convective heat fiux in 
the z direction. Substituting Eq. (10) into Eq. (8), it 
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TABLE 1 

Parameters' values of the k — e model 





Ce2 


References 


0.09 1.3 1.44 


1.92 


Pope (2000) 


1.5 


2.0 


Tennekes7l989") 



can be obtained that the temperature difference is ap- 
proximately determined by: 



d 

dxi 



dAT\ 
dxi J 



PoCp- 



dAT 

dxi 



ds 

PoToT-W, 
oz 



(11) 



where Cp is the specific heat at constant pressure. We 
omit hereafter the subscript "0" for the stationary den- 
sity and temperature, which will give rise to no confusion 
on their meanings. 

3. THE k-e MODEL FOR THE STELLAR TURBULENT 
CONVECTION 

3.1. Equations of the k-e model with buoyancy 
modifications 

The standard k-e model with buoyancy modifications 
consists of two equations (jPopd 120001 : IHossain fc Rodil 
[1982) : 

Dk d f dk\ ^ ^ 

Dt dxi (^(Te dxi) "^^^ ~'~ fc "^^ /c ' 

where k = ^UjUl is the kinetic energy of turbulence, e the 
dissipation rate of fc, and D/Dt the co-moving derivative. 
On the left hand side of Eqs. (12) and (13), the standard 
gradient diffusion hypothesis is adopted to treat the tur- 
bulent transport process, and the turbulent viscosity vt 
is approximated according to the eddy viscosity model 
as: 

vt^Ct, — . (14) 

On the right hand side of Eqs. (12) and (13), P repre- 
sents the shear production rate of turbulent kinetic en- 
ergy: 

(15, 



P = —UiUj 



while G describes the contribution from the buoyancy: 

G - (16) 

where is the temperature fluctuation of turbulence. 
There are some model parameters in above equations, 
e.g. c^, tTg, Cei, Ce2, and 0^3. Table 1 lists some choices 
of their values that are specified in applications of the 
standard k — e model. There are lots of controversy on 
the value of 0^3 , and we shall discuss this in the following 
sections. 

3.2. Steady state solutions in the fully local equilibrium 

If the left hand sides of Eqs. (12) and (13) are equal 
to zero, the turbulence is in fully local equilibrium state. 
This results in a steady state solution: 



P 

e 



Ce2 - CelCeS, 



(17) 



G 

e 



Ce2 - Cel 
Cel - CeiCeS ' 



(18) 



We use for simpli city th e turbulence parameters sug- 
gested by iTennekeai ()1989[ ) in Table 1, and define a new 
parameter as: 



Cs3 = 1 



Ce2 - Cel I 

c, 



1 



Cel 



e3 



As a result, it can be obtained that: 
^ - I ^ 



-e3 



G 

e 



C'e3 



(19) 

(20) 
(21) 



It is interesting to notice that only when 0^3 > can 
the fully local equilibrium of turbulence appear in the 
convectively unstable region. Furthermore, the condition 
that c'e3 > 1 must be satisfied to ensure a positive rate 
of shear production. 

3.3. Model for the rate of buoyancy production 

As already pointed out by iKenjeres fc Haniali^ (|2000f l 
that, in order to reproduce the ensemble roll pattern of 
convective cells, it is crucial to use an algebraic model 
for Ui-d including all its production terms: 

1\— kfTds BUi — - ^ ^\ 

1 + - Uii) = -ce- —-^ — UiUj + 4- — UjP + rjpgti)^ 
y J £ \Cp uxj oxj J 

(22) 

where the turbulent Peclet number y is defined as: 



y 



(23) 



pCp k 

On the right hand side of Eq. (22), the first term in the 
brackets is the contribution from the stratification, the 
second one from the velocity shear of convective rolling 
cells, and the last one from the buoyancy. We intro- 
duce a modification of time scale from radiative dissipa- 
tion on the left hand side (jCanuto fc Dubovikovl 119981 : 
ILi fc Yand[200l . _ 
For the auto-correlation of temperature fluctuation ■d'^ 
in Eq. (22), we adopt a local equilibrium approximation 
that assumes a complete balance bewteen its production 
and dissipation: 



-2 , kT ds 

■0 = -c e — UiV. 



(24) 



Solving for the velocity-temperature correlation wd in 
the z direction and neglecting some angle-dependent in- 
gredients, we finally approximate the rate of buoyancy 
production by the following model: 



G = — 



1 + + c'^CTT^N"^ e 
where the buoyancy frequency N is defined as: 
n T ds 



(25) 



(26) 



and the typical time scale of turbulence r is defined as: 

(27) 



_ ^ 
e ' 



It can be noticed from Eq. (25) that the property of the 
stratification determines the effect of buoyancy contribu- 
tion: G > in a convection zone where N"^ < 0, while 
G < in a stably stratified region where N'^ > 0. 

3.4. Model of the convective heat flux 

Based on similar arguments as discussed in deriving 
Eq. (25), we approximate the convective heat flux Fc 
by: 



Fc 



(28) 



with a characteristic time scale t/j that is often approx - 
imated by a model proposed by Nagano fc: KimI ()1988f ) 
as: 

(29) 



Th 



'TgT. 



In Eq. (29), the kinetic time scale of turbulence r is 
defined by Eq. (27) , and the thermal time scale of tur- 
bulence tq is approximated by: 



pCp k 
re = -: T. 

A e 



(30) 



As a result, the convective heat fiux is then expressed as: 



A 



' ,,3/2 



I3gl + y^^ + c'pCTT^N^ 



(31) 



3.5. The temperature gradient in the stellar convection 

zone 

In the stellar convective envelope, the temperature gra- 
dient V can be obtained by substituting Eq. (31) into 
Eq. (10) as: 



dlnT 
d\np 



3/2 



PgTl + y-^ + c'^,ctt'^N'^ 



7V^ (32) 



where is the radiative temperature gradient assuming 
that the heat flux is solely transferred by radiation, and 
the local pressure scale height is defined as: 



Hr, 



dr 
d\np 



P_ 

pg' 



The buoyancy frequency N can be expressed as: 

r2_ i^gT 



H 



(V-Vad), 



(33) 



(34) 



p 



where the adiabatic temperature gradient S/ad is defined 
by: 

fd\nT\ , , 

With the definition of the efficiency of convective heat 
transfer / as: 

^ V - Vad 



the buoyancy frequency can be written as: 



7V2 = 



Fin 



(V, - V„d) / - -Ej, 



(36) 



(37) 



where E' is a function of the stellar structure. 



Directly solving for / from Eq. (32) will result in a 
singular solution. Instead, we introduce a new variable: 



h = 



1-./ 



(38) 



and Eq. (32) becomes: 

(1 + y) /i^ ^ (l + 2/ - cV?/'/' - d^^CTEr'^i^ h - c' ^y''^^ = 

with Ai=c'^CTET^y + cVy^/' - y - 1 . (39) 
Equation (39) is a quadratic equation and its solution is: 



Ai + + 4 (1 + y) d^y^l-^ 

2(i + y) 



(40) 



It is easy to verify that h is always positive not only in 
the convection zone (E > 0) but also in the overshooting 
regions {E < 0). 

4. SHEAR PRODUCTION RATE FROM CONVECTIVE 
ROLLING CELLS 

4.1. Simplified structure of convective rolling cells 

The velocity structure in a convective rolling cell can 
be simplified as a circular motion in, for example, (y, z) 
plane: 

= (0, V, W) . (41) 

By use of a system of polar coordinates with the origin 
at the center of the cell, the velocities along y and z 
directions can be written as: 



V = -Vesm0 
W = Vecose, 



(42) 



where the velocity in the 6 direction can be expressed as: 

Ve = rn{r), (43) 

and the angular velocity fl is regarded as only a func- 
tion of the radius r to the center of the cell. It can be 
noticed that such a velocity structure satisfies the conti- 
nuity equation of an incompressible fluid: 



dV dW _ dn 
dy dz do 



0. 



(44) 



Based on such a simplified structure of convective 
rolling cells, the shear production rate P can be writ- 
ten as: 

ay aiy fdv dw\_ 

P = -—v^ - —w^ ^ [— + —\vw. 45) 
oy oz \oz ay J 

If the Reynolds stress is approximated by the simple eddy 
viscosity model: 



2,^ (du, , dU, 



(46) 



we obtain by substituting Eqs. (42) and (46) into Eq. 
(45) that: 

P = vtS'', (47) 

where the shear of velocity S in the convective cell is 
defined by: 
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It can be noticed that only the velocity shear contributes 
to the production of turbulence, while rigid body rotation 
has no effect on turbulence. 

4.2. The temperature difference in a convective cell 

The motion of a convective rolling cell will result in 
a temperature difference AT superimposed on the sta- 
tionary temperature background. For a circular motion 
considered above, AT approximately satisfies the follow- 
ing equation in the steady state: 



pcp 



ld_ 

r dr 



dAT\ 1 d^AT 
dr J 



Ve dAT N"^ 



de 



Ve cos 6. 

(49) 

The first term on the left hand side of Eq. (49) describes 
heat dissipation due to radiation, and the second term 
gives the temperature difference of a displacement over 
a stratified fluid. 

If the heat dissipation of radiation is ignored, the so- 
lution of Eq. (49) is: 



AT = -—rsm0 = -—z, 

I3g 



(50) 



which is exactly what the standard mixing-length theory 
adopts. The general solution of Eq. (49) can thus be 
written in the complex form as: 



AT = i—[r + B (r)] e' 



(51) 



where B is assumed to be a function of the radius r. 
Substituting Eq. (51) into Eq. (49) and introducing a 
new independent variable: 



.pCpfl 



we obtain: 



^d'^B dB , 



(52) 



(53) 



This is a first-order modified Bessel equation, and its 
solution is: 



1 

B = BohiO = Bo E u ^ni 



2m+l 



(54) 



Therefore, the solution of temperature difference can be 
approximately written as: 



AT: 



N^Rb 



r 

R'b 



Boh 



r 



(55) 



This solution satisfies the boundary condition that AT = 
at r = 0. On the other side, AT should be zero at the 
surface of a convective rolling cell, which defines a typical 
size of convective rolling cells Rb as: 



Rb 



X 



pCpQ 



(56) 



From the thermodynamic point of view, however, flu- 
ids undergoing the circular motion in a convective rolling 
cell operate as a heat engine in the convection zone, ab- 
sorbing heat when it moves inward to higher temperature 



region and releasing heat when it moves outward to lower 
temperature region and converting the difference of heat 
it absorbs and releases into its kinetic energy. As a re- 
sult, the kinetic energy of its circular motion should be 
proportional to the heat involved in it: 



2 d2 



CpAT. 



(57) 



Considering Eqs. (55), (56), and (57) and assuming the 
turbulence in the fully local equilibrium described by Eq. 

(21), we approximate the averaged angular velocity of the 
circular motion of the convective cells by: 




aHr, 



(58) 



where a is a parameter similar as in the mixing-length 
theory. 

4.3. The averaged shear of velocity in convective rolling 

cells 

Usually, the unstable stratification in the stellar con- 
vection zone is fairly weak, and the circular motion 
reaches a steady state in a convective cell if the Reynolds 
stress balances the centrifugal force of circular motion. 
With the aid of Eqs. (7) and (42), it can be obtained: 




cos^ + 



dw'^ dvw \ 
^— + I sm( 
oz 



dy J 



(59) 

Substituting the eddy viscosity model Eq. (46) into Eq. 
(59) and assuming that the shear of the circular motion 
is nearly a constant, we obtain: 



+ 



2dk 

dV 
dz 



-2 



di^dV 
dy dy 



dvt dW 

dz dz 



+ 



dW\ fdi^t .^dut . ^ 
^cos6'-h ^sm( 
dy J \dz dy 



(60) 



It is interesting to note that an axially symmetric circu- 
lar motion can induce an azimuth-dependent turbulence. 
Averaging Eq. (60) over all azimuthal angle 6 and taking 
the averaged size of convective cells into account, we can 
evaluate the averaged kinetic energy of turbulence in a 
convective cell as: 



(61) 



It can be seen that the kinetic energy is approximately 
partitioned more or less equally in the circular motion 
and in turbulence. 

The averaged shear of velocity in the convective rolling 
cells can be regarded as the ratio of their typical veloc- 
ity to their typical size. It should be recognized that 
the length scale of thermal convection is not only deter- 
mined by the averaged size of convective cells, but also 
restricted by the macro-length scale of turbulence. When 
the size of the convective cell is larger than the macro- 
length of turbulence, the fiow pattern of thermal convec- 
tion is dominated by a cellular structure fully clogged up 
with numerous convective cells. If the macro-length of 
turbulence is larger than the averaged size of convective 
cells, the structure of an individual convective cell may 
be destroyed and several convective cells will merge into 
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a larger and somewhat chaotic flow pattern with a length 
scale comparable to the macro-length of turbulence. 

Based on above arguments, we assume that the av- 
eraged shear of velocity in the convective cells can be 
approximated by: 



S-" cx 



L2' 



(62) 



where L is the macro-length of turbulence. It can be 
noticed that: 



if 
if 



Rb < i, 
Rb > L. 



(63) 



Accordingly, the velocity shear dominates the generation 
of turbulence for small convective cells, while the buoy- 
ancy production plays the major role to produce tur- 
bulence for large convective cells within them the tem- 
perature difference is able to become large enough. As 
a result, we suggest a model for the averaged shear of 
velocity in the convective cells as: 



2 q2 



1 



CL + {Ln) 



(64) 



with CL = (c^)3/4. 



4.4. Decay of the temperature difference along the 
vertical direction 

As already discussed above, the convective rolling cells 
can effectively restrict the macro-length of turbulence, 
especially for the case when the macro-length of turbu- 
lence is much longer than the size of a convective cell. It 
is therefore interesting to understand how a plume-like 
vertical flow that is formed by merging several convective 
cells restricts the development of the macro- length of tur- 
bulence in the z direction, or, in other words, how far a 
convective rolling cell can move in the vertical direction. 

The buoyancy is determined by the temperature dif- 
ference of convective flow over the stationary structure. 
We assume that the velocity of a vertical plume is char- 
acterized by: 

= (0, 0, W) , (65) 

where the velocity W in the z direction is only a function 
of distance r to axis z in the horizontal plane (x, y) as: 



W = W{x,y) = W{r), 



(66) 



which is similar to the case that axis z is not placed at 
the center of a convective cell but at the boundary of two 
adjacent convective cells. 

The equation of energy conservation Eq. (11) is now 
written as: 



A 


'1 d 




92 AT" 


pCp 


r dr 




dz'^ 



dAT 

dz 



W ^ 



I9 



The general solution of Eq. (67) is of the form: 



AT 



[z + Q (r) H (z)] 



W. 
(67) 



(68) 



where the first term in the bracket of the right hand side 
is the well known mixing-length solution that has already 
been discussed for the convective cells in Eq. (50), and 



the second term is what we are interested in at present. 
Substituting the formal solution Eq. (68) into the differ- 
ential equation Eq. (67), we obtain: 



2^ 

pCp 



Hd_ fdQ\ ^ ^d^H 
r dr \ dr J dz^ 



dz 



(69) 



Equation (69) can be solved by the method of separa- 
tion of variables: 



ld_ 

r dr 

dz^ 



1 dH 
Hb dz 



H 



= 0, 



where 



Hb 



pCpW 



(70) 
(71) 

(72) 



The solution of Eq. (70) is a zero order modified Bessel 
function, which describes an exponentially decay of ve- 
locity from the center of the plume and Rb measures 
therefore the width of the plume or the size of the con- 
vective cell. 
With the aid of a new variable: 



z = RbC, 
Eq. (71) can be written as: 



dC 



dH 

dc 



H = 0, 



where 



Rb 

""^Wb- 



(73) 



(74) 



(75) 



This is a damped oscillation equation, and its solution in 
the weak damping limit ( x ^ 1 ) is: 



H = Hne-'/^"' cos ( — 

Rb 



(76) 



It is interesting to note that the temperature difference 
along the z direction is a damped oscillating distribution, 
with the size of each unit being similar to a convective 
cell and a decaying distance of Hb- 

On the other hand, the kinetic energy of a convective 
plume is simply a result of the work done by the buoy- 
ancy: 

W"^ « -N^Rl, (77) 

because the buoyancy can only act within a distance com- 
parable to the averaged size of convective cells. Consider- 
ing Eqs. (56) and (72) and assuming that the turbulence 
is in the local equilibrium state, the decaying length of 
the temperature difference along a vertical plume can be 
approximately estimated as: 



H', 



r Hi,. 



(78) 



5. LOCAL SOLUTIONS OF THE fc-e MODEL 

5.L Choice of model parameter c'^^ 

Usually, the model parameter c^s is taken to be a con- 
stant in the literature, and choices of its value disperse 
in a wide range from -1.4 to 1.45 ([Baumert fc PetersI 
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I2000D . As already discussed in Sect. 3, the condition 
that > 1 should be satisfied to ensure a state of local 
equilibrium existing in an unstably stratified region. Ac- 
cordingly, we suggest a reasonable choice of c^3 for the 
stellar turbulent convection as: 



C'e3 = 1 



It is therefore easy to see that: 

1 < c'. < 1 + , 



e3 

CL - Cfj, 

and with the aid of Eq. (19) that: 
4 



< Ce3 < 



4 1 



3 CL 



(79) 



(80) 



(81) 



Using the model of averaged velocity shear Eq. (64) 
and substituting Eq. (79) into Eq. (20), we obtain: 



5.3. Macro-length model for stellar turbulent convection 

One of the major effects of the convective rolling cells 
is to restrict the effective length for the buoyancy do- 
ing work on the turbulent flow. Assuming that Lf, can 
be approximated by the decaying length of temperature 
difference along the z direction Hb, we obtain: 



(87) 



With the aid of Eqs. (56) and (87), we suggest a model 
for the macro-length of turbulence in the stellar convec- 
tion zone as: 



V P^P 



where a is a parameter similar as in the mixing-length 
theory. This macro-length model of turbulence will be 
applied to regions in both unstable and stable stratifica- 
tion. 



(82) 



which confirms the assumption that L in the model of av- 
eraged velocity shear is indeed the macro-length of tur- 
bulence in the state of fully local equilibrium. On the 
other hand, substituting Eq. (79) into the local equilib- 
rium condition of buoyancy Eq. (21) and using Eq. (82), 
we obtain: 

fc«-v/^V^i'^', (83) 

which is similar to the corresponding equation that has 
been used in the standard mixing-length theory. 

5.2. Macro-length model in the mixing-length theory 

In order to close the turbulent convection model dis- 
cussed above, we have to specify the macro-length of con- 
vective turbulence. From Eq. (83) we know that the pro- 
duction of turbulent kinetic energy is a direct result of 
the mechanical work done by the buoyancy. The buoy- 
ancy comes from the temperature difference between the 
turbulently convective flow and its surroundings. As a re- 
sult, the mechanical work done by the buoyancy is closely 
related to the amount of heat enclosed in the turbulent 
flow: 



c„AT « -c„-^Lb, 



(84) 



where Lf, is an effective length scale over which the buoy- 
ancy efficiently does work on the turbulent flow in the 
stellar convection zone. By use of Eqs. (83) and (84) we 
obtain: 



HriLh. 



(85) 



If the effects of the convective rolling cells are not taken 
into account, it is reasonable to assume that Lf, is just 
the macro-length of turbulence L itself. It can be seen 
from Eq. (85) that: 



(86) 



which is exactly what the standard mixing-length theory 
assumes. 



5.4. Steady state solutions in the fully local equilibrium 

In the local equilibrium state, it can be found from 
Eqs. (58), (82), and (88) that: 



V A / \c'nC'e3 



n3/2 



yO 



(89) 



Substituting Eqs. (58) and (89) into Eq. (79), the model 
parameter c'^a can be rewritten in the local equilibrium 
state as: 

c'e3 = l + T- (90) 

CL - Cf_i+y 

With the deflnition of a dimensionless quantity x as: 



Eq. (40) can be rewritten as: 




(91) 



, _ ^2 + VAi+4(l+y)c>5/2 

with A2 = c'f^crixyy^ -\- c ^y^^"^ -y - I. 

Combining Eqs. (21) and (37), the steady state solu- 
tion in the fully local equilibrium is given by: 



l + y=(cT + c',3)|y. 



x^y^ 



(93) 



By use of Eq. (32), it can be found that: 



1 + + ( 1 + 1^ j cV.//^ ^ (CT + C'.3) (94) 



Furthermore, it can be found by combining Eqs. (93) 
and (94) that: 



h=\l 



CT \ c'f.y^/'^ 



de3j l + y 



(95) 
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5.5. Comparisons with the standard mixing-length 
theory 

According to the mixing-length theory, the macro- 
length of convection is characterized by the local pres- 
sure scale height, e.g. Eq. (86). As a result, we can 
approximate the turbulent time scale r as: 



(96) 



Substituting Eq. (96) into the local equilibrium condi- 
tion Eq. (21), we obtain: 



The solution of Eq. (97) is: 



2xy = y^a;2+4(cT + c'e3)cV/ - x. 



(97) 



(98) 



On the other hand, we assume that the time scale t/,, of 
the convective heat flux is proportional to the time scale 
of turbulence r in Eq. (29). As a result, we obtain that: 



(cT + c'e3)-f = 1 

x-^y 



y 



1 



CT_ 

C'e3 



Substituting Eq. (97) into Eq. (99) we obtain: 

-1 



a:V = (cT + c',3)(l + ^) x(l-/). 

C e3 



(99) 



(100) 



Substituting Eq. (100) into Eq. (98) we finally obtain: 



yx^ + A{cT + c',3)c'^f-x^ 



1 -I- Ct/c'e3 



xil-f). 
(101) 

It can be seen that except for different model parameters, 
Eq. (101) is exactly th e same as the result of the standard 
mixing -length theory (ICox fc Giulil [19681) : 



V7" 



9 ^ 



/)■ 



(102) 



Our results are compared with that of the mixing- 
length theory in Fig. 1, in which we use h defined by 
Eq. (38) as the dependent variable in order to show more 
clearly the asymptotic behaviours of the two convection 
models. It can be seen that our model of Eq. (94) has 
the same asymptotic behaviour of a -2/3 decreasing in- 
dex as that of the mixing-length theory when a; — >■ 0, 
while declines more slowly with a -5/2 decreasing index 
than the mixing-length theory of a -4 decreasing index 
when X ^ oo. 



6. A k-. 



MODEL FOR STELLAR TURBULENT 
CONVECTION 



Based on the model assumptions discussed in the pre- 
ceding sections, the equation for the dissipation rate of 
the kinetic energy of turbulence can be written as: 



De 

m 

3 " 
2 



d 

dxi 



vt de 
Ge dxi 



P 



3 + 3(Ce3 



1) G 



e £ 
- -2— 
fc k 




Fig. 1. — The efficiency of the convective heat transfer based on 
different convection models. The solid line (MLT) represents the 
result of the standard mixing-length theory, and the dashed line 
(STC) is the result of the steady state solution of the k-e model in 
the fully local equilibrium. 

dk 

Dt dxk 



1 



2+-{c'e3-1) 



e dxk 
- d ^,c' e^L-'^k 



21^2 y^V^^^ 



e 
1 fc 



r- (103) 



There are two facts to be noticed in Eq. (103). The 
first term on the right hand side of Eq. (103) con- 
tributes to an additio nal di ff usion of turbulence as al- 
ready pointed out by iPopd ()2000f ). The second term 
provides a local equilibrium state of Eq. (103) for not 
only the unstable stratification in the stellar convection 
zone but also the stable stratification in the overshooting 
regions. In consistency to the quasi-equilibrium assump- 
tion used to derive the macro-length model of turbulence, 
we introduce a model equation for the turbulence fre- 
quency Lo as: 



Duj 



d 



Dt dx, [a^dxj VCMCe3^2 



^1/3^-4/3 



kuj'/' 



- U! 



(104) 



where the turbulence frequency u is defined by: 

_1 _ e 
T k ' 

and two quantities of the stellar structure are defined by: 



(105) 



and 



D 



H = 



pcp 
A 

aHp 



(106) 



(107) 



We adopt Eq. (90) for the model parameter c'e3 in ac- 
cordance with the quasi-equilibrium assumption used to 
derive Eq. (104). 

By use of Eqs. (37) and (40), the production rate of 
buoyancy can be written as: 



G 



2c'^Ey 



A3 + V^3 + 4^c'f,y5/^c'^CTET^y e 



(108) 
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with A3 = 1 + y + c'^y5/2 _ c^^crEr^y . 

As a result, the kinetic energy equation of turbulence can 
be written as: 



Dk d 
Dt dxi 



Vt 



dk 

dxi 



uj^ + Dfca;2 + d^,{Dkf'^uj^/'^ - d ^.crEDk 1 + Vl + ^ 



A = 



where 

Ac' f,{Dkf^^uj^/^c' f^CTEDk 



(109) 



Dkiu^ + c'f,{Dkf^^uj^/^ - d^^crEDk 



7. APPLICATIONS OF THE k-uj MODEL IN DIFFERENT 
STELLAR MODELS 

In order to compare the stellar turbulent convection 
model proposed in this paper with the standard mixing- 
length theory, we have applied both of them in calcu- 
lations of stellar evolutionary models and compared the 
results. Our stellar evolutionary models are computed 
by a stellar evolution code h04.f originally described by 
Paczy nski and Kozlowski and updated bv iSienkiewic j 
(I200I. The OPAL op acities URogers fc IglesiasI 119951: 
llglesias fc Rogerslll996[) are used in the high tempera- 
ture region and the opacities from molecules and grains 
(|Alexander &: Fergusonlfl994[ ) are used in the low tem- 
perature region. The Live rmore Laboratory equation 
of state ([Rogers et al.l 119961) are adopted in our calcu- 
lations j_Thejiuclea££ea£tioi^ upd ated accord- 
ing to IBahcall fc Pinsonneaultl (|1995D and [Harris et al.l 
(|1983f ). We have modi fied the t r eatme nt of element dif- 
fusion according to Th oul et al.l (|1994[ ). 

7.1. Comparisons of solar models 

Firstly, we apply the steady state solution (e.g. Eq. 
(94) ) of the k-uj model in the fully local equilib- 
rium to construct solar models. In order to calibrate 
our solar models to achieve the observed solar param- 
eters (solar age Tq = 4.57 Gyr, solar mass M© = 
1.9891 X 10^3 g, solar luminosity L© = 3.839 x 10^3 erg s'^ 
(|Bahcall & Pinsonneault"1995'), and solar radius i?© = 
6.9566 X 10^° cm (Habcrrcitcr & Schmutz 2008i)), we ad- 
just initial helium abundance Yq and the mixing-length 
parameter a iteratively to ensure a relative accuracy of 
10^^ for those solar parameters. The values of other pa- 
rameters of the turbulent convection model are summa- 
rized in Table 2, as well as some basic properties of our so- 
lar models. In order for our solar models to have a depth 
of the convective envelope being in agreement with the 
result of helioseismic inversions {Rcz = 0.7137?©), the 
initial metal abundance has to be fixed to Zq = 0.0208 
which results in a higher Z/ X (0.0255) compared to the 
recent observations (0.0181) (jAsplund et al.ll2009t ). 

The sound sp eed from the helioseismic inversion 
(|Basu et al.l 120091 ) is compared with those of our solar 
models with the mixing-length theory (MLT) and with 
the turbulent convection model of Eq. (94) (STC) in 
Fig. 2. It can be seen that the two solar models have al- 
most identical sound speed profiles, except for in a region 
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Fig. 2. — The sound speed difference between the hehoseismic 
inversion and solar models based on different convection models. 
The solid line (MLT) represents the result of the standard mixing- 
length theory, and the dashed line (STC) is the result of the steady 
state solution of the k-uj model in the fully local equilibrium. 
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Fig. 3. — The temperature gradient in the upper envelope of solar 
models based on different convection models. The solid line (MLT) 
represents the result of the standard mixing-length theory, and the 
dashed line (STC) is the result of the steady state solution of the k- 
uj model in the fully local equilibrium. The adiabatic temperature 
gradient is given as the dotted line for reference. 

about 5% below the solar photosphere. In order to show 
the difference of the two solar models just below the solar 
photosphere, the profiles of the temperature gradient V 
of MLT and of STC are compared in Fig. 3, along with 
the adiabatic temperature gradient Vad for reference. It 
can be seen that the turbulent convection model signif- 
icantly reduces the temperature gradient just below the 
photosphere, which is a direct result of lower decreas- 
ing rate of the convective h eat transfer efficien cy of the 
turbulent convection model ()Trampedachll20l"0l ). 

The typical velocity of turbulent motions in the con- 
vection zone is compared in Fig. 4 for the two solar 
models with MLT and STC. It can be noticed that the 
two turbulence models predict similar maximum veloci- 
ties of turbulence (about 3 x lO^cms^^) appeared just 
below the solar photosphere. But in most part of the 
convection zone, the STC predicts a typical velocity of 
turbulence (about lO'^cms"^) one order of magnitude 
lower than what the MLT does (about lO^cms"^). 
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TABLE 2 

Parameters of the convection model and properties of solar models 



Model cV CT a Yq Ys Zq Z/X Rcz(Rq) 



MLT 1.8750 0.27785 0.2474 0.0208 0.0255 0.7136 

STC 0.1 0.008 0.5 1.5536 0.27785 0.2474 0.0208 0.0255 0.7138 




igT 

Fig. 4. — The distribution of the typical velocity of turbulence 
based on different convection models in the solar convection zone. 
The solid line (MLT) represents the result of the standard mixing- 
length theory, and the dashed line (STC) is the result of the steady 
state solution of the k-uj model in the fully local equilibrium. 
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Fig. 5. — The comparison of the typical velocity of turbulence 
based on the steady state solution in the fully local equilibrium 
(STC) and on the full k-u) model including the transport effect 
with the gradient diffusion turbulence model (GDT). 



Secondly, we apply the full k-co model of Eqs. (104) and 
(109) to the whole solar interior including the convection 
zone and overshooting regions located below and above 
it, in order to study the transport effect of turbulence. 
It can be seen in Fig. 5 that the transport effect of the 
gradient diffusion turbulence model (GDT) is significant 
near the boundaries of the convection zone and in the 
overshooting regions. The typical velocity of turbulence 
roughly shows a linearly decaying law below the base 
of the solar convection zone, and according to results of 
different combinations of parameters values of the stel- 
lar turbulence model, the larger the parameter is, the 
steeper the decaying law will be, which requires a smaller 
a to calibrate the model to the solar values. However, the 
turbulent velocity extends more or less constantly into 
the solar photosphere and lower atmosphere until they 
drop abruptly to the outer boundary conditions we have 
prescribed in the upper atmosphere. It is interesting to 
note that the typical velocity of turbulence is restricted 
below Ikms^^ in the whole subphotospheric region due 
to the transport effect of turbulence, ensuring the ther- 
mally turbulent convection appearing there a subsonic 
motion as assumed in the stellar turbulence model. On 
the other hand, the typical length of turbulent motion 
remains almost a constant in the overshooting region be- 
low the base of the solar convection zone as shown in 
Fig. 6. Contrary to the local pressure scale height that 
decreases monotonically in the stellar envelope, the typ- 
ical length of turbulence increases exponentially in the 
upper convection zone and overshooting region. This is 
a direct result of radiative dissipation that effectively in- 
creases the average size of convectivc rolling cells when 
the density drops rapidly in the solar convective enve- 
lope. 
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Fig. 6.-- The comparison of the typical length of turbulence 
based on the steady state solution in the fully local equilibrium 
(STC) and on the full k-u) model including the transport eflfect 
with the gradient diffusion turbulence model (GDT). 

7.2. Comparisons of stellar evolutionary models 

We calculate a series of stellar evolution models of 

O.6M0, I.OM0, I.5M0, 3.OM0, and 7.OM0, with the ini- 
tial chemical abundance of Yf) — 0.27785 and Zq = 0.02 
and parameters of turbulence models given in Table 2. 
Element diffusion is not considered here for simplicity. 
The HR diagram is shown in Fig. 7 for the stellar mod- 
els with different masses. It can be seen that the evolu- 
tion tracks of stellar models based on the mixing-length 
theory (MLT) and on the model of stellar turbulent con- 
vection (STC) are almost identical, except that the STC 
results in a series of red giant branches with a little bit 
lower effective temperatures. 
Profiles of the typical velocity of turbulence are com- 
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Fig. 7. — The HR diagram for stellar models of different masses based on different convection models. The solid line (MLT) represents 
the result of the standard mixing-length theory, and the dashed line (STC) is the result of the steady state solution of the k-ui model in 
the fully local equilibrium. 
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Fic;. 8. — The comparison of the typical velocity of turbulence 
for an I.OMq RGB model. The solid line represents the result of 
the full k-uj model including the transport effect with the gradient 
diffusion turbulence model (GDT), and the dashed line is the result 
of the steady state solution in the fully local equilibrium (STC). 



pared for l.OM© raodels based on STC and GDT at the 
RGB bump in Fig. 8, within whieli tlie envelope eonvee- 
tion mostly penetrates inwardly into the stellar interior. 
It can be seen that the transport effect of turbulence 
is only significant below the base of the convective en- 
velope, with an almost linearly decaying law similar to 
the solar case. On the other side, the turbulent velocity 
is effectively reduced in the upper convective envelope 
and diffused as nearly a constant value of lkms~^ into 
the stellar atmosphere. It is interesting to note in Fig. 9, 
however, that the typical length of turbulence varies only 
in one order of magnitude around 10^ km in the whole 
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Fig. 9. — The comparison of the typical length of turbulence for 
an I.OMq RGB model. Others are the same as in Fig. 8. 

convective envelope, decreasing linearly below the base 
of the convective envelope and increasing rapidly in the 
stellar atmosphere. 

For a model of 7.0Mq in the main sequence, the trans- 
port effect of turbulence results in significant effect on 
the turbulent velocity in the convective core, the maxi- 
mum velocity of GDT being only 80% as seen in Fig. 10 
as that of STC. Outside the convective core, the turbu- 
lent velocity can be diffused into an overshooting region 
of about 0.5Hp wide. However, the turbulent velocity is 
decaying more and more rapidly in the overshooting re- 
gion, resulting in a considerably shorter distance of fully 
mixing for elements outside the convective core. It can 
be seen in Fig. 11, however, that the typical length of 
turbulence remains almost as a constant in the convec- 
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Fig. 10. — The comparison of the typical velocity of turbulence 
for a 7.OM0 main sequence model. The solid line represents the 
result of the full k-ui model including the transport effect with the 
gradient diffusion turbulence model (GDT), and the dashed line is 
the result of the steady state solution in the fully local equilibrium 
(STC). 
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Fig. 12. — The comparison of the typical velocity of turbulence 
for a 3.OM0 sub-giant model. The solid line represents the result of 
the full k-u) model including the transport effect with the gradient 
diffusion turbulence model (GDT) , and the dashed line is the result 
of the steady state solution in the fully local equilibrium (STC) . 
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Fig. 11. — The compajrison of the typical length of turbulence for 
a T.OMq main sequence model. Others are the same as in Fig. 10. 

tive core and in the overshooting region, and its value is 
comparable to the size of the convective core. The typical 
length of turbulence sharply increases toward the stellar 
center due to the central boundary condition, which as- 
sumes the local equilibrium of turbulence at the center 
of a star. 

Finally we apply the turbulence models of STC and 
GDT to a 3.OM0 model in the post-main-sequence evo- 
lutionary stage. According to the MLT, convection is 
underdeveloped in the envelope of such a star and tur- 
bulent motions only appear in two thin shells and have 
almost no effect on the structure of the stellar envelope. 
As seen in Fig. 12, however, the transport effect dom- 
inates the turbulent motion, and the turbulent velocity 
is diffused in the whole upper envelope upward to the 
stellar atmosphere and downward several local pressure 
scale heights below the lowest convective boundary. The 
resulted maximum velocity of turbulence is substantially 
reduced from STC's 16kms~i to GDT's 4kms-^ On 
the other hand, the typical length of turbulence remains 
almost a constant as seen in Fig. 13, and its value is 



D5 



^ 10.5 




Fig. 13. — The comparison of the typical length of turbulence for 
a 3.OM0 sub-giant model. Others axe the same as in Fig. 12. 

significantly increased to be comparable to the width of 
the convective region being formed thereabout by the 
turbulent diffusion. 

8. CONCLUSIONS AND DISCUSSIONS 

Thermal convection is characterized in stars by con- 
vective rolling cells of different scales and a fully de- 
veloped turbulence. Such semi-regular and large scale 
structures play a critical role in determining the macro 
properties of the convective motions, by both acting as 
thermal engines to transform heat into the kinetic energy 
to maintain their cellular configurations and producing 
the velocity shear within the convective rolling cells to 
generate turbulence. Due to great difficulties encoun- 
tered both physically and numerically in dealing with 
such complex flows, tm-bulent convection models, which 
are based on fully hydrodynamic moment equations and 
include many physical properties of turbulence, have not 
jet taken the contributions from the convective rolling 
cells into account. 

In order to investigate the qualitative properties of the 
convective rolling cells and to include their contributions 
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to the development of turbulence, we simplify their con- 
figurations by eddies rotating circle by circle with dif- 
ferential angular velocities, and approximately describe 
their motions by equations of mass and energy conserva- 
tion. With the aid of advanced turbulence models, this 
simple approach enable us to estimate the average size of 
the convective rolling cells by solving for the temperature 
difference over the stationary temperature background, 
and to approximate the generation rate of turbulence 
due to the average shear of velocity in the convective 
rolling cells. On the other hand, the restriction from the 
heat conversion efficiency of such convective rolling cells 
working as individual thermal engines leads to a direct 
limit on the turbulent kinetic energy enclosed in them 
and an indirect evaluation of the typical length scale of 
turbulence in the fully local equilibrium. Based on these 
considerations and model assumptions thereof, we obtain 
the static solution of the standard k-e model by properly 
introducing an important model parameter c^s, and then 
develop a k-uj model including the turbulent transport ef- 
fect to describe stellar turbulent convection. Results of 
preliminary applications to the sun and other stars with 
different masses and in different evolutionary stages are 
in good agreement with those of the standard mixing- 
length theory and numerical simulations. 

Being as an initial attempt to include the contributions 
from the convective rolling cells in the turbulent convec- 
tion models, our approach have a lot of limitations and 
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incongruities. The most significant among them lies in 
the simplistic assumption that the convective rolling cells 
are of an eddy-like configuration, which is proved to be 
inadequate for the surface layer convection in the solar 
convective envelope where the scale height of the density 
stratification is much smaller than the average size of the 
convective rolling cells. The huge variation of density in 
the convective rolling cells results in great changes on the 
topology of the flows, i.e. a distinct up/down asymme- 
try characterized by hot and slow up-flows separated by 
cold and fast down-drafts (Stein & Nordlund 1989). On 
the other hand, the properties of the convective flow is 
determined by the physical conditions along the whole 
streamline that the fluids go through during a complete 
cycle of the convective motion. Such a nonlocal nature 
of the phenomenon, along with the huge range of varia- 
tions of the density, significantly limits the validity of our 
approach by use of the local physical conditions to esti- 
mate the size of the convective rolling cells. Besides, the 
actual size and shape of the convective rolling cells may 
also be restricted by geometric conditions, for examples, 
in the stellar convective core or in thin convective shells 
of the stellar envelope with a thickness smaller than the 
local size of the convective rolling cells. 
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